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Abstract 

This the user's manual for SHARE version 2. SHARE [1] (Statistical Hadronization 
with Resonances) is a collection of programs designed for the statistical analysis of particle 
production in relativistic heavy-ion collisions. While the structure of the program remains 
similar to vl.x, v2 provides several new features such as evaluation of statistical fluctuations 
of particle yields, and a greater versatility, in particular regarding decay feed-down and in- 
put/output structure. This article describes all the new features, with emphasis on statistical 
fluctuations. 
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Program Summary 



Title of the program: SHAREv2, 



April 2006 



Computer: 

PC, Pentium III, 512MB RAM 



not hardware dependent; 



Operating system: 

Linux: RedHat 6.1, 7.2, FEDORA etc. 



not system dependent; 



Programming language: FORTRAN?? 

Size of the package: 167 KB directory, without hbraries (see 
http: / /wwwasdoc. web. cern.ch/wwwasdoc/minuit/minmain. html 
http: / / wwwasd.web.cern.ch / wwwasd / cernlib.html 
for details on library requirements ). 

Distribution format: tar gzip file 

Number of lines in distributed program (.tar hie), including test data, etc.: 26017 

Keywords: fluctuations, relativistic heavy-ion collisions, particle production, statistical models, 
decays of resonances 

Computer: Any computer with an f77 compiler 
Nature of the physical problem: 

Event-by-event fluctuations have been recognized to be the physical observable capable to con- 
strain particle production models. Therefore, consideration of event-by-event fluctuations is re- 
quired for a decisive falsification or constraining of (variants of) particle production models based 
on (grand-, micro-) canonical statistical mechanics phase space, the so called statistical hadroniza- 
tion models (SHM). 

As in the case of particle yields, to properly compare model calculations to data it is necessary 
to consistently take into account resonance decays. However, event-by-event fluctuations are 
more sensitive than particle yields to experimental acceptance issues, and a range of techniques 
needs to be implemented to extract 'physical' fluctuations from an experimental event-by-event 
measurement. 

Method of solving the problem: 

The techniques used within the SHARE suite of programs [1] are updated and extended to fluc- 
tuations. A full particle data-table, decay tree, and set of experimental feed-down coefficients are 
provided. Unhke SHAREvl.x, experimental acceptance feed-down coefficients can be entered for 
any resonance decay. 

SHAREv2 can calculate yields, fluctuations, and bulk properties of the fireball from provided 
thermal parameters; alternatively, parameters can be obtained from fits to experimental data, via 



2 



the MINUIT fitting algorithm [2] . Fits can also be analyzed for significance, parameter and data 
point sensitivity. 

Averages and fluctuations at freeze-out of both the stable particles and the hadronic reso- 
nances are set according to a statistical prescription, calculated via a series of Bessel functions, 
using CERN library programs. We also have the option of including finite particle widths of 
the resonances. A minimization algorithm, also from the CERN library programs, is used to 
perform and analyze the fit. Please see [1] for more details on these. 

Fm\iOse: 

The vast amount of high quality soft hadron production data, from experiments running at 
the SPS, RHIC, in past at the AGS, and in the near future at the LHC, offers the opportunity 
for statistical particle production model falsification. This task has turned out to be difficult 
when considering solely particle yields addressed in the context of SHAREvl.x [1]. For this reason 
physical conditions at freeze-out remain contested. Inclusion in the analysis of event-by-event 
fiuctuations appears to resolve this issue. Similarly, a thorough analysis including both fluctuations 
and average multiplicities gives a way to explore the presence and strength of interactions following 
hadronization (when hadrons form) , ending with thermal freeze-out (when all interactions cease) . 

SHAREv2 with fiuctuations will also help determine which statistical ensemble (if any), e.g., 
canonical or grand-canonical, is more physically appropriate for analyzing a given system. Together 
with resonances, fluctuations can also be used for a direct estimate of the extent the system re- 
interacts between chemical and thermal freeze-out. 

We hope and expect that SHAREv2 will contribute to decide if any of the statistical hadroniza- 
tion model variants has a genuine physical connection to hadron particle production. 

Computation time survey: 

We encounter, in the FORTRAN version computation, times up to seconds for evaluation of 
particle yields. These rise by up to a factor of 300 in the process of minimization and a further 
factor of a few when x^/NdoF profiles and contours with chemical non-equilibrium are requested. 

Accessibility: 

The program is available from: 

• The CPC program library, 

• The following website: 

http : // www . physics . arizona . edu/~ torr ieri/SHARE/share . html 

• from the authors upon request. 

SUMMARY OF NEW FEATURES (w.r.t. SHAREvl.x) 

Fluctuations: In addition to particle yields, ratios and bulk quantities SHAREv2 can calculate, 
fit and analyze statistical fiuctuations of particles and particle ratios; 

Decays: SHAREv2 has the fiexibility to account for any experimental method of allowing for 
decay feed-downs to the particle yields; 
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Charm flavor: Charmed particles have been added to the decay tree, allowing as an option study 
of statistical hadronization of J/ip, Xc, Dc, etc.; 

Quark chemistry: Chemical non-equilibrium yields for both u and d flavors, as opposed to gener- 
ically light quarks q, are considered; rj^' mixing, etc., are properly dealt with, and chemical 
non-equilibrium can be studied for each flavor separately; 

Misc: Many new commands and features have been introduced and added to the basic user 
interface. For example, it is possible to study combinations of particles and their ratios. It is also 
possible to combine all the input files into one file. 

SHARE Compatibility and Manual: 

This write-up is an update and extension of [1]. The user should consult Ref. [1] regarding the 
principles of user interface and for all particle yield related physics and program instructions, other 
than the parameter additions and minor changes described here. SHAREv2 is downward com- 
patible for the changes of the user interface, offering the user of SHAREvl a computer generated 
revised input files compatible with SHARE v2. 
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1 Introduction 



The statistical hadronization model [3-6] (SHM) assumes particles are created according to their 
phase space weight, given the locally available energy and quantum numbers. Such a reaction 
model implies that the underlying dynamics of strong interactions saturates the strength of each 
particle production quantum matrix element. 

This approach can be used to calculate the event-by-event average, as well as fluctuation 
(distribution width) and higher cumulants of any 'soft' observable. Event-by-event particle fluctu- 
ations have been subject to intense current theoretical [7-15], and experimental interest [16-19]. 
SHAREv2 will offer a standardized framework to evaluate these. 

While qualitative study of fluctuations is useful as a test of new physics, a quantitative analysis 
including both average multiplicities and their event-by-event fluctuations constitutes a powerful 
probe of hadronization conditions [20-23]. In particular, the following questions can be addressed 
when both yields and fluctuations are considered in the same model framework: 

• SHM can be falsified if and when fluctuations do not scale w.r.t. averages as expected 
in statistical physics. Moreover, only if the same set of thermal parameters gives good 
description of experimentally measured yields and fluctuations, can we claim that the SHM 
fit is physically sound. 

• As has recently been shown [24-28], the value to which the scaled variance aN (see Eq. (2)) 
for a single particle converges in the thermodynamic limit varies by as much as an order of 
magnitude when different statistical ensembles are considered. Thus, fluctuations can help 
decide if and when certain particle yields should be studied in grand canonical or canonical 
ensembles. 

• SHM fits containing both the average particle multiplicity and the fluctuation de-correlate 
the hadronization temperature T and light quark phase space occupancy 7^ (see Eq. 11 
and [1]) typical of fits when only the average multiplicities arc fitted [20]. 

Therefore, the study of both fiuctuations and yields can help to experimentally distinguish 
between the chemical equilibrium freeze-out model (T ~ 170 MeV, 7q = 1 [29]), or the best 
fit with chemical non-equilibrium at typically lower T [30]. 

• Considering the directly detectable resonance decays, fiuctuations of particle yield ratios 
offer a way to quantitatively gauge the effect of hadronic re-interactions between formation 
and thermal freeze-out [21]. 

To investigate these questions, it is necessary, in evaluation of both particle yields and fiuctu- 
ations to: 

• Incorporate all particles resonance decay trees [31] in the program structure; 

• Obtain particle yields and fiuctuations for a given set of thermodynamic parameters; 
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• i) Check if the parameters obtained by fitting particle yields are consistent with observed 

fluctuations; 

ii) Once all corrections to fluctuations due to experimental setup are understood, incorporate 
the fluctuations along with yields into the chemical freeze-out fitting procedure. 

SHAREv2 comprises a framework that addresses these challenges. 

As implied above, event-by-event particle yield fluctuations are subject to many subtle exper- 
imental effects which need to be understood and kept under control for a joint yield-fluctuation 
analysis to proceed. Further, there is the choice of statistical model ensemble in computation of 

the phase space volume: 

1) Evaluation with exact energy and discrete quantum number conservation (micro-canonical en- 
semble — MCE); 

2) In the canonical ensemble (CE), statistical energy fluctuations are allowed, conserving discrete 
quantum number (s) exactly; 

3) In the grand- canonical ensemble (GCE), statistical fluctuations of all conserved quantities occur 
— there are also mixed CE-GCE ensembles where some particle yields are conserved and other 
fluctuate. 

Clearly, the fluctuations of particle yields are most constrained in MCE and least constrained 
in GCE. Thus, although in the three ensembles, the first moments of any observable distribution, 
i.e., expectation values, coincide in the thermodynamic limit, this will not be the case for the 
fluctuations [24-26]. The choice of appropriate ensemble in the situation considered has to be 
made based on evaluation of prevailing physical conditions. 

In study of total particle yields, in the physical context of heavy ion collisions, the electrical 
charge and baryon number are fixed and, in these variables, we have to consider the CE or MCE 
if and when we are observing all particles. On the other hand, if we only observe a sub-volume 
of the system, which is exchanging energy and particles with an unobserved 'bath' consisting of 
the remainder of the reaction system, then, also conserved quantum numbers must be allowed 
to fiuctuate, which imphes use of the GCE for all observables. However, when the totahty of 
the produced particles carrying a conserved quantum number falls within the detector acceptance 
region, occasional detection of one such particle implies presence of the corresponding anti-particle, 
and thus in that case CE or MCE must be applied. 

Within the context of heavy ion physics, with reactions occurring at large energy, a study of 
fiuctuations within a narrow momentum rapidity^ acceptance window provides for the division 
between 'system' and 'bath', with the bath being the unobserved rapidity domain. In all exper- 
iments currently capable to measure fluctuations, detector acceptance is limited typically to the 
central rapidity phase space coverage. Such an acceptance domain in the boost invariant (denoted 
below as subscript b.i.) limit is equivalent to a configuration space sub-volume [32] and thus for 
both particle ratios, and particle yield width (fluctuation) [23], we have: 

^The rapidity y, defined by y = ^ In ^ , where E is the particle energy and Pz the particle momentum 

component parallel to the collision axis, is a variable additive under Lorentz transformations parallel to this axis. 
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Here, (Nj) is the event- by-event average of particle i, dNi/dy is the average number of particles 
in an element of rapidity at central rapidity. 

Similarly, barring acceptance effects (discussed in sections 1.5 and 1.6 of this work, as well as 
in [23]), the scaled variance, defined as 



_ {(AXr) _ {X^) - {Xf 
~ (x) - {X) ' 



will be given by 



da% 



(3) 

b.i. 



dy 

We conclude that, in experiments with limited central rapidity acceptance, both yields and fluc- 
tuations should be evaluated in the GCE with respect to the conserved quantum numbers (charge 
Q , baryon number 6, strangeness s — s). 

Use of the GCE, for at least some conserved charges, is required by the experimental observa- 
tion of a significant fluctuation in those charges [16-18]. This fluctuation has been found to be 
compatible with Poisson scaling, 

((AiV)2> ^ (AT) , (4) 

which is approximately followed by the GCE fluctuations. This is not the only scaling known to 
be present in this area of physics. Elementary reaction systems have been observed to follow a 
non-Poissonian scaling [33,34] w.r.t. multiplicity averages, 

{(ANr)r.{N)+c{Nf, (5) 

where c is a constant. As has been argued previously, [35-37], it is possible to describe this scaling 
by considering an extension of the Grand Canonical ensemble (variously referred to as Isobaric or 
Pressure ensemble) where system volume is also allowed to fluctuate. 

In SHAREv2, we consider only GCE yields and fluctuations and search to explore whether the 
grand canonical statistical hadronization model can quantitatively reproduce fluctuations in the 
same way as it was shown to reproduce particle yields in heavy ion A-A reactions. 



1.1 Evaluation of yields and fluctuations 

In GCE, particle yields and fluctuations can be calculated by a textbook method [38]. For a 
hadron with an energy Ep = lyp^Trnf, the energy state occupancy is, 

where the upper sign is for fermions and the lower sign is for bosons. The chemical fugacity Tj 
will be considered in section 1.2. 
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The yield average is obtained by multiplying the occupancy number Eq. 6 by the density of 
states (where V is volume and g degeneracy): 

im = gvj ^M^v)- (7) 



The fluctuation in this number is found to be: 

d{N,) 



Eqs. (7-8) can be evaluated to any desired accuracy by converting them into an expansion of 
Bessel function terms [5] , 



n=l 



ii^m - '-^t'-^i'^'rVc^), (10) 

n=l 



where W{x) — x^K2{x) (see [1], section 2, for the technical details required in doing these calcu- 
lations, as well as a discussion of particles with finite width). 



1.2 Chemical potential and chemical non-equilibrium 

The particle fugacity Tj can be obtained from the quark number content of the particle as well 
as the fugacities and phase space occupancies of the individual quark flavors. If a particle has 
Ng = N^^ + N^, Nl and A*"*, light up + down, strange and charm quarks, iV| = iVl + iVl, and 
Nl antiquarks, and isospin I^, the fugacity Tj, or, equivalently the associated particle chemical 
potential /x, will be given by: 

T, = e^^/^ = (A„7.)"^" (A.7.)''^ (A.%)^^ (Ac7c)"^^ (A.7.)''^ (A.-Tj)"'^ (A.7.)^^ (Ae-7c-)^^ , (H) 

where Xj = X^^, and 7- = 7^. The individual u, d light quark variables are related to the SHARE 
q and variables Xg — \/XuXd and A/3 = X^/Xd, see Ref. [1], and similarly 7^ = ^J%fu and 73 = 
lulldi see Eq.(23) below. The condition of chemical equilibrium, for a flavor /, imposes 7/ = 1 [38]. 
The assumption of chemical equilibrium is not automatic in a dynamically expanding system with 
a possible phase transition and, in fact, good theoretical arguments have been proposed for 7^,^ 7^ 1 
for a range of energies (see [30], and references therein). However, it is difficult, using fits to particle 
yields, to distinguish between two models based on different temperatures and 7 values. 

For instance, models based on both a higher freeze-out temperature and ^/^ = 1 [29] or a lower 
freeze-out temperature and 7^ > 1 [30] have been used to fit SPS and RHIC data. As shown 
in [23], this ambiguity is resolved when both yields and fluctuations can be considered. 
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A complication arises for hadrons such as the 7r° or the r), which are in a flavor superposition 
state. If ^u,d,s 7^ Ij the yield of the hadron with fractional flavor content is considerably altered 
by the mixing. For a meson of fractional quark number structure, 

\i >— au\uu > +ad\dd > +as\ss >, 0^ + 0^ + — 1, (12) 

the fugacity comprises the chemical yield fugacities as follows: 

= A, (7>^ + 7d«d + 7'«') ■ (13) 

Fractional flavor content has non-negligible influence on the abundances of 77° and 77', and their 
decay products, in flts which allow for chemical non-equilibrium factor 7^. The same remarks 
applies when 73 = 7m 7^^ 7^ 1 to 7r°, p°, etc. Thus, 73 7^ 1 can considerably enhance 7r° oc 73+73"^ 
yield, while maintaining tt^ yields symmetry (since 73 cancels out in tt^). 

Importantly, the evolution of quark-coalesced hadrons into flnal quark-eigenstates hadrons (like 
the oscillation of neutral kaons into K5 and K^) means that the 'source' QGP quark content will 
not, in general, be equal to the 'flnal' hadron quark content. 

1.3 Resonance decays 

Eqs. 7 and 8 can be used to calculate the event-by-event averages and fluctuations of all hadrons 
at hadronization. This, however, is quite different from the observed averages and fluctuations, 
since most hadrons are strong resonances (unstable states), which decay after freeze-out, either to 
stable particles or to other resonances. The final state particle yields can be computed by taking 
the effect of these feed-downs into account [6] . 
The ensemble average of the total yield {N^) is: 

(A^i)total = (iV,)direct + $]5,^^(iVi). (14) 

Bj^i is the probability (branching ratio) for the decay products of j to include i. 
The fluctuation after resonance feed-down is given by: 

((AiV,^0'> = B^^,{N,^, - B,^,) (TV.) + {{AN,y) . (15) 

The second term corresponds to the fluctuation in the yield of resonances. The flrst term, in the 
number of j i decays given the branching ratio bj-^i. Afj^i is the number of particles type i 
produced in the decay, so that Bj^i = Afj^i. 

Mj^i = 1 for nearly all decays of nearly all resonances; The most significant exception are decays 
to multiple 7r° s, such as 77 — > 37r°. 

1.4 Volume fluctuations and fluctuations of ratios 

The expression 10 neglects volume fiuctuations, coming from centrality cuts and dynamics of 
system expansion. These are accounted for by dividing the observed fluctuation into an extensive 
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and an intensive part, 

{{AXr) ^ {{Ax)') {Vf + {xf {{AVr) , (16) 

(x), {x') can be calculated by the statistical methods described in this section. 

It is difficult to describe the volume fluctuation coefficient {{AVY) in a model-independent 
way. The most straight forward way to deal with this problem is to choose observables insensitive 
to {{AVY). 

Any observable where {x) -C {Ax)^) would be a good candidate. This is why the fluctuation 
in electromagnetic charge has long been considered to be a promising observable [10]. 

A more general approach is to consider the event-by-event fluctuation of particle ratios [9], 
where the volume fluctuation ((AV)^) is zero by construction. Fluctuation of particles ratios 
can be calculated from the denominator's and numerator's fluctuation once the full resonance 
decay tree is known [9]. Note that, unlike in the case of particle yields, resonance decays produce 
both fluctuations and correlations, since a resonance can decay both into a numerator and a 
denominator particle. If this is the case, a high resonance admixture can considerably reduce the 
fluctuation of a ratio w.r.t. Poisson expectation. 

The formulas to be used are, for the event-by-event fluctuation of the ratio of two particles 
A^i/A^2 (for example K-/n+ [9]): 

, _ {(AN,)') , {{AN,y) JAN^AN^ 

The last correlation term is given by the resonance decay into both particles 

{AN,AN2) = (A^i A^2) - (A^i) (A^2) ^ Yl ^^'-i-^ {N)j . (18) 

j 

does not depend on the average system volume ((1^)), since it cancels between the numerator 
and the denominator. (J^i/jv2' however, does acquire a dependence on (y) since ratio fluctua- 
tions scale as {N)~^ . Hence, an analysis incorporating fluctuations of particle ratios should also 
consistently account for particle yields, and the system normalization (thermodynamic parameter 
norm, [1] section 3.1) should be considered as a flt parameter. 

The equations presented in this section can be used to compute fluctuations of particles yields 
and ratios from given SHM parameters (temperatures, chemical potentials, and phase space occu- 
pancies) . However, it has long been known [39] that fluctuations are considerably less robust than 
yields against systematic effects resulting from limited experimental acceptance. These effects, 
therefore, have to be taken into account within the SHM. The next two subsections give two such 
issues addressed within SHAREv2. 

1.5 Fluctuations and detector acceptance 

One way to separate detector acceptance effects from physics is to eliminate the former via mixed 
event techniques. A 'static' fluctuation cr^^^j is measured in a sample of fake events, constructed by 
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using tracks from different events [16]. Since tracks from different events have no correlations or 
quantum corrections, a^^^^ is determined solely by a trivial Poisson contribution as well as detector 
acceptance effects. 

Within the statistical hadronization model: 

{<J%,)stat = 1 . (19) 

For particle ratios in mixed events, the correlation term, (AA'^iAiV2) in Eq.l7, vanishes, while 
((AA^i_2)^) follow Poisson scahng. Hence, Eq. 17 reduces to 

K/iv,W = ^ + ^. (20) 

The 'dynamical fiuctuation' cr^^^ [9,39-41] corresponds to the difference between the 'raw' total 
fluctuation and the fake event fluctuation: 

^In ^ V^^ - crl^f (21) 

In certain limits, can be shown [39] to be independent of detector acceptance. Hence, com- 
paring SHM estimates of adyn to experimental measurements is more reliable than using a^. 



1.6 Correlations, resonances and detector acceptance 

Because mixed event tracks are uncorrelated, mixed event techniques cannot account for detector 
acceptance effects within particle correlations. Thus, the branching ratios appearing in Eqs. 14 
and 15 need to be modified, Bj^i — > aj^iBj^i, and in addition Eq. 17 needs to be updated. 



^m/N. - + - 2«-7]^^^y^' (22) 



here aj^i refers to the probability that particle i will end up in the detector's acceptance region 
given particle j is in the region, while ai2 measures the probability that both decay products will 
be inside this region. 

For a boost invariant azimuthally complete system, aj^i — 1 since particles leaving the detector 
acceptance region will be balanced by particles entering. However, in general ai2 < 1, since if 
a resonance is outside the detector acceptance region both particles can not be inside it, and the 
intrinsic particle decay momentum adds a rapidity scale to the system, breaking boost invariance 
[20]. 

See [20] for an illustration of how to calculate 0:12. While such a comprehensive calculation 
is outside the scope of the current version of the program, we offer the user the possibility of 
entering an 0:12 for any resonance decay as an input parameter. Sec section 3 on how to do this. 
In practice, this should only be necessary for a few most frequent and energetic resonance decays, 
such as p — >• TTTT and K* — >• Ktt. 
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2 Implementation of GCE fluctuations in SHAREv2 



Experimental event-by-event fluctuation data points were implemented in the SHARE interface 
in a similar manner as yield and ratio data points (see [1], section 3.4). The tag which denotes 
that a fluctuation is being calculated is fluct_yld. A statement such as 

particlel fiuct_yld data Agtat Asyst fit? 
will calculate aj^ of particlel (defined in Eq. 2). 

If fit? is set to 1, this data point is used within a fit together with the experimentally measured 
value data and the statistical (Agtat) and systematic (Agyst) error . The format of the data- line 
is exactly the same as in SHAREvl.x ( [1] section 3.4). 

To calculate the fluctuation of a ratio, particlel should be substituted by the data point 
number where the ratio is defined. For instance, if the 5th data point (from the top) is a K~ /t:~ 
ratio, than it's fluctuation is given by: 

05 fiuct_yld data A^tat A^yst fit? 

(See section 5.2 of this paper for a more general treatment of this data point referencing). 

SHAREv2 implements most definitions of dynamical fluctuations used to date by experimental 
collaborations. These are implemented as additional tags of fiuct_xxx type, where xxx refers to 
different ways the experimental measurement is presented. 

The possible types of data points are: 

fiuct_dyn To calculate Udyn — -x/cr^ — ^Itat^ ^ rneasured in [41], 
fluct_dnr To calculate adyn = — — as suggested in [9]. 

'^stat 

3 Decay feed-down and particle yields 

3.1 Particle decay acceptance data files 

As shown in section 1, decay feed-down is a fundamental component of the statistical hadronization 
model. However, the limited coverage of most detectors means that the feed-down coefficients will 
acquire an experimental correction, corresponding to the probability that the decay products of 
a given resonance formed within the detector acceptance region will also be in that region. As 
shown in section 1.6, these corrections need to be considered when calculating both fluctuations 
and yields. 

Weak decays, such as A — pn' (most protons at RHIC are, in fact, given by feed-down from 
hyperons) , are particularly susceptible to experimental acceptance, as they occur at a macroscopic 
distance from the primary vertex. Hence, weak experimental feed-down corrections include a 
geometrical as well as a momentum space component. 

Since the 'parent' particles are not always directly observed, SHARE must be able to compute 
final hadron multiplicities including experimental feed-down coefficients for all decays where this 
effect is non-negligible. 
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SHAREvl.x allowed the user to input experimental (weak) feed-down contributions to pro- 
duced particle yields via four acceptance coefficients: 

Ks anything, Kl anything, Y Mesons, Y — > baryons (see [1], section 3.4.1). 

It turns out this approach was not sufficiently flexible: for instance, undetected E — > A^tt 
feed-down can be treated very differently from A — > jm" corrections, considering the difference in 
lifespan, and (vertex) acceptance cuts applied. 

Moreover, the experimental acceptance of different hyperon nucleon weak decays, such as 
S~ — s> Att^ as compared to A ^ p-K^ is likely to be considerably different. Finally, different weak 
decays of the same hadron can have varying acceptances, compare Kl — > 37r, with Kl — > Trez/, and 
with Kl — > 7r//i/. A similar acceptance problem may arise in special cases involving strong decay 
chains when the acceptance region is particularly narrow and/or the particles rapidity distribution 
is not well understood. 

A more ffexible way of treating weak decay contributions to particle yields is therefore nec- 
essary. Speciflcally, there should be an easy way to allow for any arbitrary decay/reaction con- 
tributing to any data point. SHAREv2 provides such a possibihty through user defined decay 
feed-down files. 

In data file containing the experimental results to be fitted (see [1] section 3.4.3), a weak decay 
control file is now signaled by a statement of the type: 
Weakdecay File. feed 

where File. feed is a 9-letter filename. The program then obtains the decay acceptance weights 
from File. feed, an ASCII file in a format similar to the decay tree files (described in [1] section 3). 

Fig. 1, and the attached input files provided with the SHARE package, show how to implement 
the weak decay acceptance coefficients. While many weak feed-down files might be involved in 
the same analysis, generally, they are experiment-specific, and hence can be kept track of in 
a systematic way. Alternatively, all weak feed-down files and experimental data-files can be 
combined in a single large file, using the methods described in section 5.1 of this paper. 

In more detail, a typical line in a feed-down file will be: 
Parent Daughteri Daughter2 all/lst/2nd/cor coeff 
or, for 3-body decays, 

Pcirent Daughteri Daughter2 Daughters all/lst/2nd/3rd/cor coeff 
The switch all/lst/2nd/3rd refers to the daughter to which the decay coefficient applies. 

all means that the decay coefficient is the same for all daughters, while lst/2nd/3rd means 
only the lst/2nd/3rd daughter will be removed from the experimental yield. For example, in the 
A — > 7rp decay in STAR [42,43], STAR accepts the nucleon from the A decay but not the tt, and 
this fine tuning of the decay is clearly quite important as a relatively large fraction of all nucleons 
comes from weak A decays. 

cor refers to the fractional contribution of the acceptance to the two particle correlation 
{Daughteri Daughter2) induced by a common resonance decay from parent, denoted as in 
Eq.22 (section 1.6). 

SHARE will renormalize the decay PcU-ent — > all/lst/2nd/3rd/cor 
by the coefficient coeff when calculating all data points after the given weeikdecay statement. 

It is therefore possible to assign a different weak decay file to each data point, or assume that a 
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Example of Experimental data file: 



#STAR data-points 

# 

weakdecay star. feed 

# 

Ka0492plu pi0139plu 0.156 0.0208 0. 

Ka0492min pi0139min 0.15 0.02 0. 1 

#PHENIX data-points 

# 

weakdecay phen.feed 

# 



Example of weak feed-down file: 



# 

#K_L-> pi 
Ka04921ng 
Ka04921ng 
Ka04 921ng 
Ka04 921ng 
#K_L-> 3pi 
Ka04921ng 

# 

Ka0492sht 

# 

Lml 1 1 5 zer 

#Lmlll5zer 

# 

Lml 1 1 5 z rb 
#Lmlll5zrb 
# 



-K L corrections 



lepton 

pi 1 3 9min 
pi0139plu 
pi 1 3 9min 
pi0139plu 

pi0139plu 



pi0139plu 



elOOOOplu 
elOOOOmin 

muOOOOplu 
muOOOOmin 



nue zer 
nueOOOzer 

numO zer 
numO zer 



all 
all 

all 
all 



pr0938plu 
neO 93 9zer 

pr0938plb 
ne0939zrb 



pi0139min pi0135zer all 
K_S corrections 

pi0139min all 0. 
Hyperon corrections 

pi0139min 2nd 0. 
pi0135zer 2nd 0. 



pi0139plu 
pi0135zer 



2nd . 
2nd 



Figure 1: An example of the SHAREv2 weak feed-down acceptance coefficient implementation. 

group of data points are subject to the same set of weak decay yield contribution, (e.g., in general 
all data points from the same experiment should have the same weak decay file). The way to do 
this is the same as in vl.x ( [1] section 3.4.3): when the program reads a weeikdecay statement, 
it assigns the current decay pattern to each data point encountered until a new wecikdecay 
statement is met. 

Two special case exist, for which no File. feed file is needed: 
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weakdecay UNCORRECT 'uncorrected' (from the perspective of experimental data set) means 
that aU weak decays contributions to particle yields are fully accepted by SHAREv2. 
weakdecay NOWK_FEED means that all particle yields are computed without contributions 
from weak decays, from the perspective of experimental data this means that either all weak 
decay products are not accepted and/or have been all corrected for in experimental yields, as, 
e.g., apphes to some NA49 results. 

When fluctuations arc considered, it is important to deal carefully with experimental correc- 
tions, which arc neither close to total (close 100 % accepted) nor null (close to % accepted). Weak 
decay corrections of the daughter particles are usually correlated with each other in momentum 
space, so the straight-forward application of Eq. 8 will not be a good description of fluctuations 
with a non-trivial detector acceptance function [13]. In this case, it is better to use dynamical, 
rather than total, fluctuations as discussed in section 1.5. 

3.2 Compatibility with SHAREvl.x experimental data files 

The improved weak decay treatment docs not impair compatibility of experimental input files be- 
tween SHAREv2 and SHAREvl.x. SHAREv2 will read a SHAREvl.x experimental data-file, and 
automatically calculate applicable contributions for each weak decay based on the information con- 
tained in the SHAREvl.x weakdecay statement. A line will be printed within the sharerun.out 
output file that signals a SHAREvl.x format weeikdecay statement was encountered. 

In addition, an output v2 weakdecay file called we£Lk#vl.x (where # refers to the data point 
number) is automatically generated translating the vl.x weak decay information into v2 format. 
The user is advised to eventually change all weakdecay lines to 
weakdecay weak^^vl.x 

as the v2 format is considerably more powerful and less amenable to systematic error stemming 
from an incomplete understanding of weak decays. 

4 Quark chemistry 

SHAREvl.x input files listed particle chemical content by total isospin / and its third component, 
/s, as well as the number of light {q, either u or d), strange (s) and charm (c) quarks (See section 
3.2 of [1]). In SHAREv2, u and d quarks are now separately accounted for. The particle listing 
format is: 

name iiiciss width spin I 13 u d s au ad as c ac MC 
where name is the particle's 9-character name, / and are the total and third component of the 
isospin, u, d, s, c are the numbers of up, down, strange and charm quark numbers, while au, 
ad, as, ac are the respective anti quark numbers. The format of the table is otherwise identical 
to that discussed in [1], section 3.2. 

To check for the possibility that phase space occupancy differs for the up and down quarks, a 
statistical model fit parameter (see [1], section 3.1) gamS (73) has been introduced, such that: 

lu = 7973, Id = 7g/73- (23) 
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The quark/anti quark numbers can be fractional, to account for the superposition states described 
at the end of section 1.2. 

To calculate {u, d, s) quark abundance in the statistically hadronizing QGP system (different, 
in general, from the freeze-out content, as shown in section 1.2), new bulk variables tot_u_qgp, 
tot_d_qgp and tot_s_qgp were introduced. These can be used in the same way as other bulk 
variables (see [1], section 3.4.2). 

4.1 Charm mesons 

Charmed particles have now been added in the files peirticles.data and partnowdt.data. Their 
nomenclature follows the general structure as described in [1], section 3.2. Dcxxxxxxx refer to 
Dc mesons, Dsxxxxxxx, chixxxxcc to Xc states and psixxxxcc to J/ip states. 

Their abundance is regulated by the chemical potential Ac and the phase space occupancy 7c, 
described in [1], section 3.1. 

5 User Interface files and new commands 
5.1 New single file control 

As described in detail in [1], section 3, SHARE relies on quite a few input files. 

• The run-file shcirerun.data 

• The particles fist 

• The particles decay tree 

• The initial values of the thermodynamic parameters 

• The experimental data points 

• Initialization for each fit parameter 

This structure makes it easy to quickly explore regions of parameter space within an analysis in 
progress. However, this system makes it easy to mistakenly lose a successfully completed and 
saved analysis, since a change in each of the files could considerably alter the end result. The 
introduction of weak decay correction file (see section 3) aggravates this problem. 

SHAREv2, therefore, makes it possible to combine some, or all, input files into a single file. 
Once the user found an optimum analysis, all input files involved in it can be combined into one 
large sharerun.data file, which can be easily kept for future reproduction and modification. 

This is done by changing the extension (.data or .feed) of the filename into .HERE. If the 
program encounters a filename ending in .HERE, it assumes the relevant input is immediately 
following the given line within the currently read file. The subsequent format is assumed to be 
unchanged from what it would have been had a separate file been opened (comments, etc.). The 
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#Many-files format 
# 

READ THERMJNI th_neq.data 



#One-file format 
# 

READ THERM_INI th_neq.HERE 
# — content of thermal file starts 
temp 0.14 



READ TOTALDATA tot200mix.data 



<within the tot200mix.data file> 
weakdecay star.feed 



accu 0.01 

# — content of thermal file ends 
READ TOTALDATA tot200mix.HERE 
#— content of experiment file starts 
pi0139plu prt_yield 286.4 24.2 0. 1 

weakdecay star.HERE 

# — content of weak decay file starts 



Ka492sht pi0139plu pi0139min all 0.7 



# — content of weak decay file ends 



# 

CALC FITRATIOS fitnw20M._neq 



# — content of experimental file ends 

# 

CALC FITRATIOS fitnw20M._neq 



Figure 2: Left: sharerun.data calling other input files. Right: One-file format. 

only difference is that a * symbol on a new line has to be present at the point where the separate 
file would have ended. When the program encounters the * symbol, it switches back to reading 
the 'earlier' file, that is prior to the insert .HERE. 

See Fig 2, and the provided file shcirerun.data_onefile, for an example of how this works. 

5.2 Combining data points 

SHAREv2 gives the possibility to refer to a different data point within the given fit, and/or 
combine two data points, in order to fit the sum or a product of two particles. This feature was 
described, in the case of fluctuations of ratios, in section 2 of this paper. 

The referring data point consists of one or two (for a combination) numbers, corresponding to 
the position, in the input file, of the point (s) being referred to. 
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Two numbers united by an operation sign (+, — , X, /) will add, subtract, multiply and divide 
two data points. For instance, if the first data points from the top of the file (see [1] section 3.4, 
for a detailed explanation of the format) are: 



LmlllSzer 


pi0139min 


Data 


^stat 


^syst 


Fit? 


LmlllSzrb 


pi0139plu 


Data 


^stat 


^syst 


Fit? 


then. 












01X02 


prt_yield 


Data 


^stat 


^syst 


Fit? 



wiU fit (AA)/(7r+7r-), while 

01 fiuct_dyn Data Agtat \yst Fit? 



will fit the dynamical A/tt fluctuation, as described in section 2. 

To flt (A + A)/(7r+ + TT") (but NOT the separate yields), the input flle wiU read: 



LmlllSzer 


prt_yield 


Data 


^stat 


^syst 





LmlllSzrb 


prt_yield 


Data 


Astat 


A syst 





pi0139plu 


prt_yield 


Data 


^stat 


Agyst 





pi0139min 


prt_yield 


Data 


^stat 


^syst 





01+02 


03+04 


Data 


^stat 


^syst 


1 



NOTE: SHARE was written in FORTRAN??. Feature mentioned in this subsection use implicitly 
recursive code. SHAREv2 has been tested on several compilers and platforms, and found to 
work. However, compilers and operating systems vary — we would like to know if and when you 
experience problems. 

5.3 Miscellaneous 

The following (small) modifications were made in SHAREv2 compared to SHAREvl: 

5.3.1 Expanded parameter set 

The expanded parameter set includes as noted before, Eq. (23), gam3 which allows to incorporate 
a different u, (i-flavor phase space occupancy. A further new variable dvol describes statistical 
pressure ensemble fluctuations in volume (Section 1, Eq. (16)). The provided input flle sets and 
flxes dvol to zero and gam3 to unity, since experimental measurements sensitive to these param- 
eters have not as yet been published. 

All details about how to configure these parameters, and fix or relax them in the context of 
fits to experimental data, are unchanged w.r.t. vl.x, described in [1] sections 3.1, 3.6 and 3.? 

5.3.2 Data point sensitivity analysis 

Command DFIT, within the file sharerun.data can turn on and off the given data point as a 

point to be fitted. 

The syntax for this command is 

DFIT [Datapoint n.] [Fit(0/1)] 

where Datapoint n. refers to the data point's position in the experimental data file from the 
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top, while Fit(0/1) turns this point on (1) or off (0) as a point to be fitted. For instance, the 
following input in sharerun.data: 
READ TOTALDATA tot200mix.data 
DFIT 5 1 

CALC FITRATIOS fltnw20M-kpi 
DFIT 5 

CALC FITRATIOS fitnw20M-nkpi 

performs two fits. 

The first, saved in file fitnw20M-kpi uses the 5th data point in tot200mix.data when calculating 
the (to be minimized). 

The second one, saved in file fitnw20M-nkpi, does not. 

5.3.3 Data point sensitivity profiles 

Command SNSPROFIL calculates the data point sensitivity. The sensitivity is defined as the ra- 
tio between the data point's SHM prediction for a given statistical parameter, and SHM prediction 
at the best fit value for that parameter. 

The syntax of SNSPROFIL is the same as DATPROFIL in [1], section 4. The two com- 
mands operate in the same way: all parameters, except the one on the abscissa, are minimized at 
each point in the profile. 
Thus, the command 

CALC SNSPROFIL temp 0.1 0.2 100 5 

will calculate a sensitivity profile for the temperature, going from 0.1 to 0.2 GeV, with 100 points, 
of the fifth data point within the experimental data-file. 

5.3.4 Additional output in and statistical significance profiles 

profile commands now output the following files: 

name. log A fit output for each point in the profile, in the same format as the usual fit output 
file ( [1], Fig. 4). 

name.chi2, name.stsg Commands SNSPROFIL and DATPROFIL also output the pro- 
file (extension *.chi2) and Ptme profile (extension *.stsg). 

5.3.5 Improved treatment of fit errors 

SHAREv2 automatically runs the MINOS algorithm [2] if the fit can not get a robust estimate of 
the errors. This results in a considerable improvement of error-handling. This update entails no 
changes in the user interface or output format. 
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6 Comparison with previous versions 



6.1 Testing SHAREv2 

SHAREv2 was extensively tested for programming and physics errors: 

• SHM Calculations and fit results for SPS and RHIC energies were verified to be equal 
between SHAREv2 and SHAREvl.x reference results; 

• SHAREv2 reads SHAREvl.x weak decay input. The equivalence between the two treat- 
ments, when weak decay files are designed to reproduce SHAREvl.x format, was shown to 
all decimal places; 

• Fluctuations of conserved quantities (such as ((AQ)^)) were compared before and after 
resonance decays. The conservation of this quantity implies that the enhanced fiuctuations 
after all resonances decayed are exactly balanced by multiplicity correlations between the 
resonance decay products. This holds true to two decimal places (up to two-step correlations 
arising from decays such as ^(1600) K{892)n — >• Kin:. These correlations are not tracked 
by SHARE, but their contribution is below 1%). 

6.2 SHAREvl.x bugs found 

While developing and testing SHARE v2, several minor bugs and choice issues were found in the 
previous version SHAREvl.x, The most noteworthy issues which lead to sometimes noticeable 
(beyond line width) changes in the results are: 

SHAREvl.x, vl.2 The Bessel function series was incorrectly truncated for large ■jg (close to 
pion B-E condensation); 

SHAREvl.3 Quark flavor mixing error in calculation of mesons such as r] and 0, for jq ^s] 

SHAREvl.l-vl.3 The most relevant issue is actually not an error but lack of versatihty in the 

handling of E — > pvr decays: E-particles decay weakly, like the As and the Ss. However, 
unlike A and S, S-decays are not experimentally reconstructible since at least one of S-decay 
products is neutral. In general, SHAREvl.x particles from these decays were included in 
the yield count. However, it turned out that while some experiments had much less than 
full acceptance for these decays, other experiments, e.g., NA49, have removed E — > p feed- 
down via Monte-Carlo simulations accounting for the experimental acceptance of the decay 
products, with E yields obtained from the observed A-yields. 

Working with patched SHAREvl.x. we realized that E-decay issue mattered in that some fits 
got better allowing for a modified E weak decay pattern. Issues, such as this one, prompted us to 
introduce a more general treatment of weak decay particle yield contributions in SHARE v2. 
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7 Installation 

The^ SHAREv.2.1 program code and input files are contained in a tar.gz archive. The file 
sharev2.1.tar.gz is available at 

http : //www . physics . arizona . edu/~ torrier i/ SHARE/ share . html 

To unpack it, create a SHARE directory, put the archive in it, and execute the following com- 
mands: 

gunzip sharev2.1.tar 
tcir -xf shcirev2.1.tcir 

The following files will then be created, enough for a complete 'representative' run of SHARE: 

decays. data The complete Particle Data Group decay tree (section 3.3 in [1]); 

dec_no.data An empty decays file, useful for testing the program calculations (abundancies re- 
duce to modified Bessel functions) as well as studying the role of resonances in stable particle 

ratios; 

fortrat A shell script compiling (in f77) the FORTRAN code which should be modified depending 
on location of FORTRAN (g77 or f77 or f95) and CERN library of programs; 

particles. data Particle properties, with full widths (section 3.2, in [1])); 

partnowdt.data Particle properties, with no widths. Calculations with this input file require 
considerably less computational time, and it suffices when there are no resonances in the fit; 

ratioset.data The FORTRAN fit input file (section 3.5, in [1]); 

sharerun.data A 'representative' run input file (section 4, in [1]); including an analysis of fiuc- 
tuations and yields similar to what was presented by members of our collaboration in [21]; 

shcLrerun.data_onefile The same file in the single file format, as explained in section 5.1; 

samplefit200 A directory containing the output files generated by running the provided 'share- 
run.data', as a debugging/comparison standard; 

sharev2.1.f SHAREvl.l FORTRAN code. The header contains information about bug fixes; 

thermo.data A representative thermal parameter input file (section 3.1, in [1]). It is set to 
reasonable non-equilibrium fit values; 

tot bcLr200. data A representative data input file (section 3.4, in [1]) containing ratios, yields and 
fiuctuations drawn from RHIC experiments, as of July 2004 (see references in [21]); 

steir.feed An example of a decay feed-down coefficients file. See section 3. 

^This chapter is nearly identical to the READMEv2.html file found on the SHARE webpage 
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Note that SHARE requires CERN libraries, to be downloaded separately from 
http: / /ww w as d. w eh. cern.ch/ www as d/ cernlih.html 

The compiler statement (in file 'fortrat') is 
f77 -L/usr/local/cern/pro/lib -o sharev2.1.exe sharev2.1.f -Imathlih -Ikernlih -Ipacklib -C 
which assumes that the CERN hbraries are in directory 
/usr /local /cern/pro /lib 

If this is not true on your system, fortrat should be changed accordingly. 

Once the directory is unpacked, the program should be compiled with 
./fortrat 

After this, typing 
. /shcire v2 . 1 .exe 

should produce a correct run with a detailed output which shows the program's capability. A copy 
of these files, produced on our computer system, is also included in sharev2.1.tar.gz within a 
directory called samplefit200 for comparison with the files produced by the installed program. 
Several output files are produced, with the following names as default. The contents of each file 
are explained in detail in [1] section 4. 

fit*. out Fit output files; 

graph* Fit output graphics (experiment, fitted values, calculated values); 
prof* profiles and correlation functions for the various fits. 
See [1] about more details about these files's contents. 

8 Status, conclusions, and future plans 

One of the areas of current intense interest in the field of high energy heavy ion reactions is the 
understanding of the mechanisms of soft hadron production (chemical freeze-out), that is the study 
of how the energy confined in the central fireball turns into matter in a multi particle production 
process. 

The SHARE suite of programs is an analysis tool of particle yields addressing the following 
questions: 

• What is the chemical freeze-out temperature, chemical potentials and volume? 

• What are the physical properties of the fireball which hadronizes? 

• Is the hadron system in chemical equilibrium at freeze-out? 

• What is the degree of re-interaction between hadronization and freeze-out? 

The need for SHARE arises from recognition that the book keeping task involved in the correct 
application of the statistical hadronization model is considerable, often transcending the resources 
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available to individual researchers. The current SHARE v2 program follows on SHAREvl.x [1] 
adding three significant novel features: a) flexible handling of particle decay feed-down, b) fluctu- 
ations, and c) complete u, d, s, c flavor content treatment. We note that, since SHAREvl.x was 
released, another analysis package appeared, THERMUS [44]. THERMUS has more convenient 
setup for the experimental data analysis environment. On the other hand, THERMUS is not 
addressing many of the features SHARE offers, including chemical non-equilibrium and now in 
the current SHAREv2, fluctuations, and light flavor details. 

SHARE allows an analysis of experimental data that can address indirectly questions related to 
the dynamics of the flreball evolution, since chemical non-equilibrium implies a fast hadronization. 
However, the analysis of particle yields cannot give results of greater precision than is inherent 
in the data it treats, and this is not yet good enough to resolve the question about the nature of 
SHM as seen in the recent references [29,30,45-47], and the possible relation to the onset of a 
phase transition, [48]. SHAREv2 offers the additional analysis feature, the fluctuation in hadron 
yields and ratios, which will help to settle these issues, should precise particle yield experimental 
data not become available. Moreover, SHAREv2 through the study of the consistency of particle 
yield and fluctuation, can test the SHM in depth. 

The development of phenomenological tools capable of falsifying statistical hadronization mod- 
els is, of course, far from over. Possible extensions of chemical freeze-out model, in future version 
of SHARE, might include a canonical ensemble module, allowing to test SHM chemical non- 
equilibrium in small physical systems, and the introduction of an opacity parameter to correct the 
yields of observed resonances [49]. Another possible future development would entail extending 
SHARE towards a detailed description of momentum distributions. This would be somewhat dif- 
ferent from THERMINATOR [50] model which relies on scaling in rapidity, and thus, only applies 
to ultra high energy collisions. We note that particle momentum distributions are dependent, in 
addition to the physics incorporated in SHAREv2 also on the dynamical evolution of the emitting 
source, and on the degree of resonance rescattering following on the chemical freeze-out. Strategies 
how these two interwoven effects can be disentangled are being developed. For further discussion 
of these and related questions we refer to a recent review [51]. 
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